path_proj = here::here()
path_source = file.path(path_proj, "source")

# source
source(file.path(path_source, "simulation", "simulations_functions_final.R"))
## 
## Attaching package: 'dplyr'
## The following objects are masked from 'package:stats':
## 
##     filter, lag
## The following objects are masked from 'package:base':
## 
##     intersect, setdiff, setequal, union
## 
## Attaching package: 'tidyr'
## The following object is masked from 'package:reshape2':
## 
##     smiths
## The following object is masked from 'package:rstan':
## 
##     extract
# source from other paper
source(file.path(path_source, "funcs_german_paper", "plotReportingTriangle.R"))
# data import
SARI <- as.data.frame(read.csv(file.path(path_proj, "data", "raw", "clean_data_srag_epiweek_delay_table_PR.csv"),))

# set week
set_week_start("Sunday")
SARI$date <- get_date(week = SARI$epiweek, year = SARI$epiyear)
fixed_q <- file.path(path_proj, "source", "models",
                     "trunc", "1.stan_model_fixed_q_trunc.stan")
fixed_b <- file.path(path_proj, "source", "models",
                     "trunc", "2.stan_model_fixed_b_trunc.stan")
b_poly <- file.path(path_proj, "source", "models", "trunc",
                     "3.stan_model_b_poly_trunc.stan")
b_spline <- file.path(path_proj, "source", "models", "trunc",
                     "4.stan_model_b_spline_trunc.stan")

compiled_models <- list(
  fixed_q = cmdstan_model(fixed_q),
  fixed_b = cmdstan_model(fixed_b),
  b_poly = cmdstan_model(b_poly),
  b_spline = cmdstan_model(b_spline)
)
source(file.path(path_source, "functions", "plot_function.R"))
source(file.path(path_source, "functions", "fit_function.R"))

# setting: delay,  number of days
seed <- 123
D <- 15
N_obs <- 200

# try-on in region 41002
data_41002 <- SARI %>% filter(regionalsaude == 41002) %>%
  select(-c("epiweek","epiyear","regionalsaude")) %>%
  relocate(Notifications, .after = last_col())
rownames(data_41002) <- data_41002$date

#transfer to cumu matrix
data_41002[,c(1:27)] <- cumulative_matrix(as.matrix(data_41002[,c(1:27)]))

# now <- as.Date("2024-02-01")
scoreRange <- seq(as.Date("2009-07-20"),as.Date("2009-08-30"),by="7 day")

# input the true case vector
case_true <- as.matrix(data_41002[, c("Notifications")])
rownames(case_true) <- as.character(data_41002$date)

out <- nowcasting_moving_window(data_41002[, c(1:21)], scoreRange = scoreRange, 
                                case_true = case_true,
                                start_date = as.Date("2009-07-01"),
                                #start_date = scoreRange[1] - weeks(20),
                                D = 20, sigma_b = 0.5, seeds = 123,
                                models_to_run = c("fixed_q", "fixed_b", "b_poly","b_spline"),
                                compiled_models = compiled_models,
                                iter_sampling = 5000, iter_warmup = 1000, refresh = 0,
                                num_chains = 3, suppress_output = T)
## ====================
## now=2009-07-20 (1/6)
## ====================
## Running MCMC with 3 sequential chains...
## 
## Chain 1 finished in 0.8 seconds.
## Chain 2 finished in 0.9 seconds.
## Chain 3 finished in 0.9 seconds.
## 
## All 3 chains finished successfully.
## Mean chain execution time: 0.9 seconds.
## Total execution time: 2.9 seconds.
## 
## Running MCMC with 3 sequential chains...
## 
## Chain 1 finished in 0.3 seconds.
## Chain 2 finished in 0.3 seconds.
## Chain 3 finished in 0.3 seconds.
## 
## All 3 chains finished successfully.
## Mean chain execution time: 0.3 seconds.
## Total execution time: 1.2 seconds.
## 
## Running MCMC with 3 sequential chains...
## 
## Chain 1 finished in 0.6 seconds.
## Chain 2 finished in 0.3 seconds.
## Chain 3 finished in 0.4 seconds.
## 
## All 3 chains finished successfully.
## Mean chain execution time: 0.4 seconds.
## Total execution time: 1.6 seconds.
## 
## Running MCMC with 3 sequential chains...
## 
## Chain 1 finished in 0.4 seconds.
## Chain 2 finished in 0.5 seconds.
## Chain 3 finished in 0.5 seconds.
## 
## All 3 chains finished successfully.
## Mean chain execution time: 0.5 seconds.
## Total execution time: 1.7 seconds.
## 
## ====================
## now=2009-07-27 (2/6)
## ====================
## Running MCMC with 3 sequential chains...
## 
## Chain 1 finished in 1.4 seconds.
## Chain 2 finished in 1.2 seconds.
## Chain 3 finished in 1.1 seconds.
## 
## All 3 chains finished successfully.
## Mean chain execution time: 1.2 seconds.
## Total execution time: 4.0 seconds.
## 
## Running MCMC with 3 sequential chains...
## 
## Chain 1 finished in 0.4 seconds.
## Chain 2 finished in 0.4 seconds.
## Chain 3 finished in 0.6 seconds.
## 
## All 3 chains finished successfully.
## Mean chain execution time: 0.4 seconds.
## Total execution time: 1.6 seconds.
## 
## Running MCMC with 3 sequential chains...
## 
## Chain 1 finished in 0.5 seconds.
## Chain 2 finished in 0.5 seconds.
## Chain 3 finished in 2.3 seconds.
## 
## All 3 chains finished successfully.
## Mean chain execution time: 1.1 seconds.
## Total execution time: 3.5 seconds.
## 
## Running MCMC with 3 sequential chains...
## 
## Chain 1 finished in 0.6 seconds.
## Chain 2 finished in 0.7 seconds.
## Chain 3 finished in 0.6 seconds.
## 
## All 3 chains finished successfully.
## Mean chain execution time: 0.6 seconds.
## Total execution time: 2.1 seconds.
## 
## ====================
## now=2009-08-03 (3/6)
## ====================
## Running MCMC with 3 sequential chains...
## 
## Chain 1 finished in 2.9 seconds.
## Chain 2 finished in 2.5 seconds.
## Chain 3 finished in 2.2 seconds.
## 
## All 3 chains finished successfully.
## Mean chain execution time: 2.5 seconds.
## Total execution time: 7.8 seconds.
## 
## Running MCMC with 3 sequential chains...
## 
## Chain 1 finished in 1.7 seconds.
## Chain 2 finished in 1.6 seconds.
## Chain 3 finished in 1.4 seconds.
## 
## All 3 chains finished successfully.
## Mean chain execution time: 1.6 seconds.
## Total execution time: 5.0 seconds.
## 
## Running MCMC with 3 sequential chains...
## 
## Chain 1 finished in 1.9 seconds.
## Chain 2 finished in 2.2 seconds.
## Chain 3 finished in 1.9 seconds.
## 
## All 3 chains finished successfully.
## Mean chain execution time: 2.0 seconds.
## Total execution time: 6.4 seconds.
## 
## Running MCMC with 3 sequential chains...
## 
## Chain 1 finished in 2.2 seconds.
## Chain 2 finished in 3.0 seconds.
## Chain 3 finished in 2.5 seconds.
## 
## All 3 chains finished successfully.
## Mean chain execution time: 2.6 seconds.
## Total execution time: 8.1 seconds.
## 
## ====================
## now=2009-08-10 (4/6)
## ====================
## Running MCMC with 3 sequential chains...
## 
## Chain 1 finished in 3.3 seconds.
## Chain 2 finished in 3.3 seconds.
## Chain 3 finished in 3.1 seconds.
## 
## All 3 chains finished successfully.
## Mean chain execution time: 3.2 seconds.
## Total execution time: 9.9 seconds.
## 
## Running MCMC with 3 sequential chains...
## 
## Chain 1 finished in 3.2 seconds.
## Chain 2 finished in 2.6 seconds.
## Chain 3 finished in 2.6 seconds.
## 
## All 3 chains finished successfully.
## Mean chain execution time: 2.8 seconds.
## Total execution time: 8.9 seconds.
## 
## Running MCMC with 3 sequential chains...
## 
## Chain 1 finished in 5.2 seconds.
## Chain 2 finished in 4.6 seconds.
## Chain 3 finished in 5.0 seconds.
## 
## All 3 chains finished successfully.
## Mean chain execution time: 5.0 seconds.
## Total execution time: 15.3 seconds.
## 
## Running MCMC with 3 sequential chains...
## 
## Chain 1 finished in 4.6 seconds.
## Chain 2 finished in 4.4 seconds.
## Chain 3 finished in 4.4 seconds.
## 
## All 3 chains finished successfully.
## Mean chain execution time: 4.5 seconds.
## Total execution time: 13.8 seconds.
## 
## ====================
## now=2009-08-17 (5/6)
## ====================
## Running MCMC with 3 sequential chains...
## 
## Chain 1 finished in 4.1 seconds.
## Chain 2 finished in 4.2 seconds.
## Chain 3 finished in 4.5 seconds.
## 
## All 3 chains finished successfully.
## Mean chain execution time: 4.3 seconds.
## Total execution time: 13.1 seconds.
## 
## Running MCMC with 3 sequential chains...
## 
## Chain 1 finished in 3.7 seconds.
## Chain 2 finished in 4.4 seconds.
## Chain 3 finished in 4.2 seconds.
## 
## All 3 chains finished successfully.
## Mean chain execution time: 4.1 seconds.
## Total execution time: 12.6 seconds.
## 
## Running MCMC with 3 sequential chains...
## 
## Chain 1 finished in 9.5 seconds.
## Chain 2 finished in 9.4 seconds.
## Chain 3 finished in 10.3 seconds.
## 
## All 3 chains finished successfully.
## Mean chain execution time: 9.7 seconds.
## Total execution time: 29.6 seconds.
## 
## Running MCMC with 3 sequential chains...
## 
## Chain 1 finished in 5.3 seconds.
## Chain 2 finished in 6.0 seconds.
## Chain 3 finished in 7.8 seconds.
## 
## All 3 chains finished successfully.
## Mean chain execution time: 6.4 seconds.
## Total execution time: 19.3 seconds.
## 
## ====================
## now=2009-08-24 (6/6)
## ====================
## Running MCMC with 3 sequential chains...
## 
## Chain 1 finished in 4.5 seconds.
## Chain 2 finished in 4.5 seconds.
## Chain 3 finished in 4.4 seconds.
## 
## All 3 chains finished successfully.
## Mean chain execution time: 4.5 seconds.
## Total execution time: 13.8 seconds.
## 
## Running MCMC with 3 sequential chains...
## 
## Chain 1 finished in 5.7 seconds.
## Chain 2 finished in 5.2 seconds.
## Chain 3 finished in 5.8 seconds.
## 
## All 3 chains finished successfully.
## Mean chain execution time: 5.6 seconds.
## Total execution time: 16.9 seconds.
## 
## Running MCMC with 3 sequential chains...
## 
## Chain 1 finished in 15.2 seconds.
## Chain 2 finished in 13.9 seconds.
## Chain 3 finished in 14.9 seconds.
## 
## All 3 chains finished successfully.
## Mean chain execution time: 14.6 seconds.
## Total execution time: 44.2 seconds.
## 
## Running MCMC with 3 sequential chains...
## 
## Chain 1 finished in 7.5 seconds.
## Chain 2 finished in 7.4 seconds.
## Chain 3 finished in 7.7 seconds.
## 
## All 3 chains finished successfully.
## Mean chain execution time: 7.5 seconds.
## Total execution time: 22.9 seconds.
results <- nowcasts_plot(out, models_to_run = c("fixed_q", "fixed_b", "b_poly", "b_spline"))

results$nowcasts
## [[1]]
##         date case_true case_reported mean_fixed_q lower_fixed_q upper_fixed_q
## 1 2009-07-05        30            17    133.31839      46.64999      401.1522
## 2 2009-07-12        64            11     73.65665      22.21361      230.0231
## 3 2009-07-19       223             2     56.81406      11.07223      192.6796
##   mean_fixed_b lower_fixed_b upper_fixed_b mean_b_poly lower_b_poly
## 1    19.802307     11.818367      42.16699    17.78680    11.230400
## 2    10.504008      4.755791      25.90012    14.33943     4.949041
## 3     7.518189      1.797386      20.99066    15.02397     2.345150
##   upper_b_poly mean_b_spline lower_b_spline upper_b_spline
## 1     33.31173     17.182216      11.998743       24.95523
## 2     38.86748      8.730459       4.407380       15.99467
## 3     43.36155      5.577643       1.322419       14.54677
## 
## [[2]]
##         date case_true case_reported mean_fixed_q lower_fixed_q upper_fixed_q
## 1 2009-07-05        30            18    111.82049      45.33370      274.6593
## 2 2009-07-12        64            11     65.84152      24.38255      167.1970
## 3 2009-07-19       223            18     91.26846      32.52952      234.8767
## 4 2009-07-26       639            32    312.52320     100.94102      845.6314
##   mean_fixed_b lower_fixed_b upper_fixed_b mean_b_poly lower_b_poly
## 1     19.98174     13.893788      31.25885    19.41739    12.739805
## 2     11.42998      6.725125      19.56273    10.07463     6.006593
## 3     15.43489      8.418845      28.31064    12.32245     7.201909
## 4     46.22332     23.913853      94.03930    31.72390    18.992755
##   upper_b_poly mean_b_spline lower_b_spline upper_b_spline
## 1     43.40384      19.15419      14.016965       27.46154
## 2     20.93934      11.27083       6.629580       19.90693
## 3     25.21783      15.02600       8.223486       27.25052
## 4     69.53310      41.05433      24.780068       69.41412
## 
## [[3]]
##         date case_true case_reported mean_fixed_q lower_fixed_q upper_fixed_q
## 1 2009-07-05        30            19    106.16023      51.06794      231.3775
## 2 2009-07-12        64            25     89.31174      40.93426      193.8198
## 3 2009-07-19       223            91    336.06359     159.10863      732.6172
## 4 2009-07-26       639           329   2885.77450    1310.22925     6407.2708
## 5 2009-08-02      1331           272  17519.04340    7537.87350    39700.5975
##   mean_fixed_b lower_fixed_b upper_fixed_b mean_b_poly lower_b_poly
## 1     312.5958      64.33369      1275.000   16487.280    1229.1370
## 2     272.9074      53.62035      1132.745    5555.680     461.5512
## 3    1015.7924     198.15353      4235.884    7117.365     672.1575
## 4    6533.0157    1269.40375     26963.682   16105.715    1621.8465
## 5   14554.9245    2734.95650     61123.365   13406.400    1292.1498
##   upper_b_poly mean_b_spline lower_b_spline upper_b_spline
## 1     78793.81      22.62987       16.05667        36.6624
## 2     26405.45      43.88084       16.64506       122.0198
## 3     33576.62     295.06132      113.89650       709.6286
## 4     74073.96    1251.80798      625.86873      2495.7258
## 5     61660.34     721.40568      374.07145      1678.9457
## 
## [[4]]
##         date case_true case_reported mean_fixed_q lower_fixed_q upper_fixed_q
## 1 2009-07-05        30            22     101.1988      52.48704      206.0062
## 2 2009-07-12        64            39     115.6864      59.28950      235.2852
## 3 2009-07-19       223           154     532.3109     276.01545     1087.5330
## 4 2009-07-26       639           463    3296.0519    1691.39250     6701.4555
## 5 2009-08-02      1331           563    7883.8010    4012.39450    16188.9775
## 6 2009-08-09      1324           183    6701.1398    3329.85925    13867.4200
##   mean_fixed_b lower_fixed_b upper_fixed_b mean_b_poly lower_b_poly
## 1     664.5866      139.2546      2720.393   20933.383    1743.4123
## 2     772.1524      161.4177      3151.671    9917.254     941.5613
## 3    3460.7446      726.0558     13839.747   17488.672    1864.3603
## 4   17821.7100     3702.5860     72474.202   36744.624    4038.6457
## 5   35965.4325     7432.5652    146923.800   31145.257    3272.3892
## 6   23525.8658     4791.7868     95778.742    8879.344     835.0984
##   upper_b_poly mean_b_spline lower_b_spline upper_b_spline
## 1    112129.92      22.80778       17.03424       34.29437
## 2     52240.59     113.12357       34.57167      354.13070
## 3     91383.51     976.56980      385.76797     2345.92350
## 4    193812.32    3290.25938     1550.31900     6773.90100
## 5    157226.37    1877.17243      968.19653     4044.52550
## 6     45372.17     779.78248      269.35700     2563.51500
## 
## [[5]]
##         date case_true case_reported mean_fixed_q lower_fixed_q upper_fixed_q
## 1 2009-07-05        30            23     90.37367      49.29757      177.1377
## 2 2009-07-12        64            47    123.33489      66.48288      245.2817
## 3 2009-07-19       223           179    564.83845     311.26555     1107.5337
## 4 2009-07-26       639           519   2711.89706    1482.69625     5357.3300
## 5 2009-08-02      1331           749   5461.32794    2986.19050    10787.9850
## 6 2009-08-09      1324           549   5180.85360    2823.86550    10193.5200
## 7 2009-08-16       860           145   4031.76048    2121.54850     8135.5607
##   mean_fixed_b lower_fixed_b upper_fixed_b mean_b_poly lower_b_poly
## 1     825.9353      144.0060      3325.223   21696.877    1252.5477
## 2    1152.8408      197.1500      4616.323   13952.930    1047.3467
## 3    5277.4571      907.5011     20748.145   29849.471    2828.3058
## 4   23824.1118     4061.5960     93118.340   65002.842    6927.0462
## 5   46710.7421     7938.0132    184770.125   62972.679    6733.1358
## 6   43012.6666     7214.7577    171637.475   29594.031    2971.3160
## 7   25582.5461     4163.9907    102242.300    9233.209     800.6328
##   upper_b_poly mean_b_spline lower_b_spline upper_b_spline
## 1    123191.45      23.07897       17.81705       31.51302
## 2     80784.31     238.92735       61.40564      776.20422
## 3    163114.75    1689.78263      653.41512     3887.71025
## 4    352419.95    3467.46140     1824.79325     6798.57050
## 5    345831.37    2994.99293     1557.71750     6429.31850
## 6    163362.07    2762.66170     1295.12625     5732.20225
## 7     51634.50     918.25508      230.27995     3504.31475
## 
## [[6]]
##         date case_true case_reported mean_fixed_q lower_fixed_q upper_fixed_q
## 1 2009-07-05        30            25     81.99376      48.45588      148.7370
## 2 2009-07-12        64            53    125.33546      73.32550      224.3693
## 3 2009-07-19       223           198    559.78898     329.34952     1004.1915
## 4 2009-07-26       639           550   2327.84895    1367.77950     4189.2705
## 5 2009-08-02      1331           932   4500.98974    2638.69725     8079.0722
## 6 2009-08-09      1324           767   4384.64140    2566.93750     7900.2925
## 7 2009-08-16       860           417   3373.63354    1967.41375     6129.6857
## 8 2009-08-23       879           249   5833.38616    3325.97500    10632.5375
##   mean_fixed_b lower_fixed_b upper_fixed_b mean_b_poly lower_b_poly
## 1     377.9449      104.0189      1499.828    12431.61     617.1515
## 2     590.4376      160.0075      2288.200    10063.84     694.3955
## 3    2667.9433      722.2312     10318.495    24016.91    2204.7238
## 4   10961.1402     2948.7370     42443.912    53547.40    6062.0782
## 5   21754.0638     5807.2428     84881.142    59098.46    7723.4870
## 6   21515.0723     5683.8283     84707.110    33273.22    4601.1780
## 7   16059.0498     4158.0995     63611.677    14523.50    1896.2160
## 8   21205.7208     5344.4838     84359.415    11451.78    1311.5397
##   upper_b_poly mean_b_spline lower_b_spline upper_b_spline
## 1     67606.31      26.62506       18.72072       31.99615
## 2     54624.59     341.40807       93.16127      956.24750
## 3    126258.27    1605.09631      746.39595     3145.93750
## 4    266306.37    2292.88952     1303.80900     4530.89975
## 5    287455.47    3114.13039     2019.61825     5240.50500
## 6    165668.87    3835.33043     2008.81750     7637.41075
## 7     74280.96    2552.98415     1036.18700     5705.89825
## 8     58840.54    1297.60192      390.42020     4189.35175
results$plots
## [[1]]

## 
## [[2]]

## 
## [[3]]

## 
## [[4]]

## 
## [[5]]

## 
## [[6]]

nowcasts_plot(out, models_to_run = c("fixed_q"))[["plots"]]
## [[1]]

## 
## [[2]]

## 
## [[3]]

## 
## [[4]]

## 
## [[5]]

## 
## [[6]]

nowcasts_plot(out, models_to_run = c("fixed_b"))[["plots"]]
## [[1]]

## 
## [[2]]

## 
## [[3]]

## 
## [[4]]

## 
## [[5]]

## 
## [[6]]

nowcasts_plot(out, models_to_run = c("b_poly"))[["plots"]]
## [[1]]

## 
## [[2]]

## 
## [[3]]

## 
## [[4]]

## 
## [[5]]

## 
## [[6]]

nowcasts_plot(out, models_to_run = c("b_spline"))[["plots"]]
## [[1]]

## 
## [[2]]

## 
## [[3]]

## 
## [[4]]

## 
## [[5]]

## 
## [[6]]